/*
 * Open Source Physics software is free software as described near the bottom of this code file.
 *
 * For additional information and documentation on Open Source Physics please see:
 * <http://www.opensourcephysics.org/>
 */

/**
 * particles start in a left box and
 * eventually fill the simulation cell showing the approach to equilibrium. There
 * is a button which allows you to reverse the velocities and see the particles go
 * back to their original positions.
 * @author Jan Tobochnik
 * @author Peter Sibley
 */
package org.opensourcephysics.stp.neto;

import org.opensourcephysics.display.Drawable;
import org.opensourcephysics.display.DrawingPanel;

import java.awt.*;

public class LJgas2boxBeta1 implements Drawable {
    public double x[], y[], vx[], vy[], ax[], ay[], alpha[];
    public Color color[];
    public int numberOfParticles = 270;
    public double cellLength; //tamanho da caixa
    public double rho = 0.2; // fixed density
    public double x0;
    public double temperature = 5.0;
    public int steps = 0;
    public double dt = 0.01;
    public int n0 = 0, n1 = 0;
    double radius = 0.3;
    public double cut2 = 3.0 * 3.0;
    public double t = 0;
    public double vAngular = 10; //velocidade angular gerada pelo B
    public double vRadius=4; //raio do movimento circular de B

    /**
     * Constructor Diamagnetism
     */
    public LJgas2boxBeta1() {
        setArrays();
    }

    /**
     * Metodo principal que faz os calculos a cada passo
     * */
    public void step() { // Velocity Verlet algorithm
        double dt2half = 0.5 * dt * dt;
        double dthalf = 0.5 * dt;
        n0 = 0;
        n1 = 0;
        for (int i = 0; i < numberOfParticles; i++) {
            x[i] += vx[i] * dt + ax[i] * dt2half;
            y[i] += vy[i] * dt + ay[i] * dt2half;
            x[i] = image(x[i]);
            y[i] = image(y[i]);
            vx[i] += ax[i] * dthalf;
            vy[i] += ay[i] * dthalf;
            if (x[i] < x0) { //Faz a contagem para saber se esta do lado direito ou esquerdo
                n0++;
            } else {
                n1++;
            }
        }
        calcRotation();
        wallColision();


        /*accel(); //calcula a mudan�a na aceleracao devido a colisao entre as particulas
        for(int i = 0; i<numberOfParticles; i++) {
         vx[i] += ax[i]*dthalf;
         vy[i] += ay[i]*dthalf;
       }             */
        steps++;
        t += dt;      
    }

    /**
     * Este metodo calcula a acelera��o caso ocorra uma colis�o entre as particulas.
     * Nos demais casos, a acelera��o � zero.
     * */
  public void accel() {
    double dx, dy, fx, fy, r2, fOverR, rm2, rm6;
    for(int i = 0; i<numberOfParticles; i++) {
      ax[i] = 0;
      ay[i] = 0;
    }
    for(int i = 0; i<numberOfParticles-1; i++) {
      for(int j = i+1; j<numberOfParticles; j++) {
        dx = pbc(x[i]-x[j]);
        dy = pbc(y[i]-y[j]);
        r2 = dx*dx+dy*dy; //raio ao quadrado
        //esse trecho cuida da colis�o  
        if(r2<cut2) {
          //System.out.println("colisao");
          rm2 = 1.0/r2;
          rm6 = rm2*rm2*rm2;
          fOverR = 48.0*rm6*(rm6-0.5)*rm2;
          fx = fOverR*dx;
          fy = fOverR*dy;
          ax[i] += fx;
          ay[i] += fy;
          ax[j] -= fx;
          ay[j] -= fy;
          //PEsum += 4*(rm6*rm6 - rm6);
          //Vsum += dx*fx + dy*fy;     //virial calculation
        }
      }
    }
  }

    /**
     * Isso deve ser periodic boundary conditions, utilizado para determinar a colisao
     */
    private double pbc(double s) {
        if (s > 0.5 * cellLength) {
            s -= cellLength;
        } else if (s < -0.5 * cellLength) {
            s += cellLength;
        }
        return s;
    }


    /**
     * Esse metodo que cria a parede infinita
     */
    private double image(double s) {
        if (s > cellLength) {
            s -= cellLength;
        } else if (s < 0) {
            s += cellLength;
        }
        return s;
    }

    private void wallColision() {
        for (int i = 0; i < numberOfParticles; i++) {
            //verifica se o x atingiu alguma das paredes
            if (x[i] + radius > cellLength || x[i] - radius < 0) {
                vx[i] = -vx[i];
                alpha[i] += Math.PI;
            }
            //verifica se o y atingiu alguma das paredes
            if (y[i] + radius > cellLength || y[i] - radius < 0) {
                vy[i] = -vy[i];
                alpha[i] += Math.PI;
            }
        }
    }

    private void calcRotation() {
        for (int i = 0; i < numberOfParticles; i++) {
            alpha[i] += vAngular * dt;
            //System.out.println("x[i]" + vRadius * dt * Math.cos(alpha[i]));
            //System.out.println("y[i]" + vRadius * dt * Math.sin(alpha[i]));
            x[i] += vRadius * dt * Math.cos(alpha[i]);
            y[i] += vRadius * dt * Math.sin(alpha[i]);
        }
    }

    public void reverse() {
        for (int i = 0; i < numberOfParticles; i++) {
            vx[i] = -vx[i];
            vy[i] = -vy[i];
        }
    }

    public void draw(DrawingPanel myWorld, Graphics g) {
        if (x == null) {
            return;
        }
        int pxRadius = Math.abs(myWorld.xToPix(radius) - myWorld.xToPix(0));
        int pyRadius = Math.abs(myWorld.yToPix(radius) - myWorld.yToPix(0));
        if (pxRadius < 1) {
            pxRadius = 1;
        }
        if (pyRadius < 1) {
            pyRadius = 1;
        }
        for (int i = 0; i < numberOfParticles; i++) {
            int xpix = myWorld.xToPix(x[i]) - pxRadius;
            int ypix = myWorld.yToPix(y[i]) - pyRadius;
            g.setColor(color[i]);
//      Color color = Color.BLACK;
//      if(x[i]>cellLength/2) color=Color.BLUE;
//      g.setColor(color);
            g.fillOval(xpix, ypix, 2 * pxRadius, 2 * pyRadius);
        }
        g.setColor(Color.black);
        int xpix = myWorld.xToPix(0);
        int ypix = myWorld.yToPix(0);
        int lx = myWorld.xToPix(cellLength) - myWorld.xToPix(0);
        int ly = myWorld.yToPix(cellLength) - myWorld.yToPix(0);
        g.drawRect(xpix, ypix, lx, ly);
//    g.drawLine(xpix+lx/2, ypix, xpix+lx/2, ypix+ly); linha que divide no meio
        g.drawLine(xpix, ypix, xpix, ypix + ly);
        g.drawLine(xpix + lx, ypix, xpix + lx, ypix + ly);
        g.drawLine(xpix, ypix, xpix + lx, ypix);
        g.drawLine(xpix, ypix + ly, xpix + lx, ypix + ly);
    }


    // ---------------- Inicializacao ----------------------------------------------------------------------------------
    
    public void setArrays() {
        int N = numberOfParticles;
        x = new double[N];
        y = new double[N];
        vx = new double[N];
        vy = new double[N];
        ax = new double[N];
        ay = new double[N];
        alpha = new double[N];
        color = new Color[N];
    }

    public void initialize() {
        cellLength = Math.sqrt(numberOfParticles / rho);
        System.out.println("cellLength=" + cellLength);
        reset();
        setArrays();
        setPositions();
        setVelocities();
        setAngles();
        //accel();
    }

    public void setPositions() {
        boolean overlap;
        //double ds = cellLength/2.0; isso faz preencher apenas metade da caixa
        double ds = cellLength; //preenche a caixa inteira
        for (int i = 0; i < numberOfParticles; ++i) {
            do {
                overlap = false;
                do{
                    x[i] = ds * Math.random();
                } while (x[i] < radius || x[i] > cellLength - radius); //garante que a bola nao fica presa na borda

                do{
                    y[i] = cellLength * Math.random();
                } while (y[i] < radius || y[i] > cellLength - radius); //garante que a bola nao fica presa na borda

                int j = 0;
                while ((j < i) && !overlap) { //isso evita que as bolas j� comecem colidindo
                    double dx = x[i] - x[j];
                    double dy = y[i] - y[j];
                    dy = pbc(dy);
                    if (dx * dx + dy * dy < 1.0) {
                        overlap = true;
                    }
                    j++;
                }
            } while (overlap);
        }
    }

    private void setAngles() {
        for (int i = 0; i < numberOfParticles; ++i) {
            alpha[i] = Math.random() * 2.0 * Math.PI;
            color[i] = new Color((float) Math.random() ,(float) Math.random(),(float) Math.random());
        }
    }

    public void setVelocities() {
        double tpi = 2.0 * Math.PI;
        for (int i = 0; i < numberOfParticles; ++i) {
            double r1 = Math.random();
            double r2 = Math.random() * tpi;
            vx[i] = Math.sqrt(-2 * temperature * Math.log(r1)) * Math.cos(r2); // sets random initial velocities
            vy[i] = Math.sqrt(-2 * temperature * Math.log(r1)) * Math.sin(r2);
        }
        // zero Momentum
        double vxsum = 0.0;
        double vysum = 0.0;
        for (int i = 0; i < numberOfParticles; ++i) {
            vxsum += vx[i];
            vysum += vy[i];
        }
        double vxcm = vxsum / numberOfParticles;
        double vycm = vysum / numberOfParticles;
        for (int i = 0; i < numberOfParticles; ++i) {
            vx[i] -= vxcm;
            vy[i] -= vycm;
        }
    }

    public void reset() {
        steps = 0;
        t = 0;
        n0 = 0;
        n1 = 0;
        x0 = cellLength / 2.0;
    }

    public void zeroAverages() {
        steps = 0;
        n0 = 0;
        n1 = 0;
    }
}

/*
 * Open Source Physics software is free software; you can redistribute
 * it and/or modify it under the terms of the GNU General Public License (GPL) as
 * published by the Free Software Foundation; either version 2 of the License,
 * or(at your option) any later version.

 * Code that uses any portion of the code in the org.opensourcephysics package
 * or any subpackage (subdirectory) of this package must must also be be released
 * under the GNU GPL license.
 *
 * This software is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston MA 02111-1307 USA
 * or view the license online at http://www.gnu.org/copyleft/gpl.html
 *
 * Copyright (c) 2007  The Open Source Physics project
 *                     http://www.opensourcephysics.org
 */
